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We discuss particle production mechanisms for heavy ion collisions. We present an 
argument demonstrating how the fluctuations of the number of produced particles in a series 
of classical emissions can account for KNO scaling. We predict rapidity correlations in the 
particle production in the event by event analysis of heavy ion collisions on the rapidity 
scales of the order of l/a s . 

I. INTRODUCTION 

In ultrarclativistic heavy ion collisions in the upcoming experiments at the Relativistic Heavy Ion 
Collider (RHIC) at Brookhaven and Large Hadron Collider (LHC) at CERN, a very high multiplicity of 
particles will be produced. Most of the models predict the multiplicities of produced hadrons per unit of 
rapidity in a given event to be of the order of Q| 

dN o 

where the exact number depends on the particular model. For these high multiplicities the scale of statistical 
fluctuations in one unit of rapidity is expected to be of the order of a ~ W dN/ dy so that 

Thus the multiplicities would not change much within each narrow bin in rapidity, as long as the number 
of particles produced in that rapidity bin is large. Collisions may therefore be studied on an event by event 
basis with little ambiguity. The fluctuations of multiplicities in each of the rapidity bins could be observed 
by examining several different events. One might wonder whether the numbers of particles produced in 
different rapidity bins are completely independent of each other, or there are certain correlations among 
them. In this paper we are going to address the question of rapidity correlations in the framework of 
McLerran- Venugopalan model of particle production by a series of classical emissions |^| . 

In McLerran- Venugopalan model pj the high multiplicity per unit area gives rise to a dimensionfull 
parameter characterizing the nuclear collision 
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with R the nuclear radius. It was argued that this scale represents the typical transverse momentum of the 
particles produced in the early stages of the nuclear collisions . Due to the high multiplicities of Eq. (|l|) 
this typical transverse momentum scale of produced partons given by Eq. (j^) may become large compared 
to the QCD scale 

A 2 » A 2 QCD - (4) 

The QCD coupling at the scale A is weak, a s (A) -C 1. Together with the fact that the number of color 
charge sources in the nuclei is large this allows us to assume that the gluons produced in a nuclear collision 
could be described by the classical field of the colliding nuclei Here by classical field we mean that 

the gluon field is a solution of Yang-Mills equations of motion with the nuclei providing the source term 
for the equations |^-|4|]. A renormalization group (RG) approach has been developed recently to account 
for a series of classical emissions H|). The procedure invented in ||[| is the following: at each step of the 
evolution in rapidity each nucleus through classical emissions produces several partons with soft longitudinal 
momenta. At the next step of the RG the partons produced in the previous step are included in the source, 
off which the next generation of even softer partons is emitted. That way the produced soft partons modify 
the density of color sources, for which one writes a renormalization group equation |.|. The exact for m of 
this equation is not important for our discussion here. We just note that at the lowest order it reproduces 
the well known Balitsky, Fadin, Kuraev, Lipatov (BFKL) equation || and the full equation should be able 
to resum multiple reggeon exchanges in the structure functions. The crucial assumption which we have to 
keep in mind is that the series of classical emissions in both nuclei provide us with a set of color charge 
sources, which give rise to the color field of the produced gluons. This is illustrated in Fig. [j] The QCD 
evolution, which consists of real terms (classical emissions) and virtual corrections provide us with the color 
sources, which are denoted by crosses in Fig. |l|. 
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FIG. 1. An example of generation of color charge sources through the QCD evolution in a nucleus. The multiple 
reggeon evolution, which is described in the text, after a number of iterations, that also include emissions of gluons, 
produces a number of sea gluons in the rapidity region far away from the nucleus, which later on act as sources of 
color charge for classical gluon field. The sources are denoted by crosses. 

To construct the classical field one has to solve Yang-Mills equations treating the color charges generated 
through evolution of Fig. [l] in both nuclei as contributions to the source term in the equation. For the field 
to be purely classical the generated color sources should be separated by the rapidity interval Sy < l/a s . 
An example of the diagram contributing to the classical field generated this way is given in Fig. |^. The 
classical field of colliding nuclei corresponds to gluon production in the approximation, where one resums all 
powers of the parameter a 2 s L ||, where L is the number of sources at a given impact parameter. If the total 
rapidity interval between the colliding nuclei is not very large, Y ~ l/a s , so that the quantum evolution has 
not become important yet, the valence quarks in the nucleons will be the sources of color charge and the 
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resummation parameter will become a 2 s A 1 / 3 j|. One would have the nucleons instead of the crosses in Fig. 
^. Finding the classical field even for this somewhat simplified situation appears to be a complicated task. 
The problem has been solved analytically only for the gluon production in the case of a proton scattering 
on a nucleus in Several attempts have been made for the case of gluon production in nucleus-nucleus 
collisions Q|, giving only the lowest order in a s result. There have also been done extensive numerical studies 
of the nucleus-nucleus collisions [ [To| . In this paper we will not be interested in the detailed structure of 
the field. What will be important to us is the fact that the classical field is boost invariant, and, therefore, 
is rapidity independent. Thus if one fixes the configuration of sources in Fig. || the classical field in the 
rapidity interval Ay ~ l/a s would give rise to a rapidity independent distribution of the produced particles 
within this rapidity interval. 

The classical gluon field has another interesting property. In the regime when the number of color 
sources is large, making the resummation parameter seizable, a 2 L ~ 1, the classical field of Fig. ^ becomes 
strong, ~ X/g ||,||,f7|,|]] . Since gluons are massless bosons and the typical phase space density is large, we 
expect that the gluon distribution function at momentum scales less than A will saturate with occupation 
number §-@,fHJ 

1 dN 

l/a s . (5) 



irR 2 cPprdy 



That way the field is very strong and non-perturbative, even though it could be obtained in the weak coupling 
regime by usual perturbative methods. The corresponding multiplicities are large [Eq. (0)]. 
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FIG. 2. A classical gluon field produced in a nuclear collision as envisaged in the text. Through QCD evolution 
nucleons in nuclei give rise to color charge sources (partons) denoted by crosses (see Fig. [j]). The sources interact 
with each other to produce the gluon field. 

If the field is classical in each collision, the produced particles arise from a fixed source on an event 
by event basis. Each event is characterized by a certain configuration of the color charges in the source, 
off which the classical field of Fig. || is emitted. The source results from the fluctuations in the color 
charge density of those quarks and gluons at rapidities larger or smaller than that at which we measure the 
distribution of produced particles. We know that the spectrum of fluctuations is Poissonian for a coherent 
state corresponding to a fixed source density which produces the field. This is similar to the results of 
reggeon calculus: after fixing the sources the classical fields might be considered as reggeons independently 
producing particles (see Fig. ||). The spectrum therefore will be Poissonian |l6|-|26[]. Moreover, the classical 
effective action in a random background field of describes both the dynamics of the gluon fields and the 
fluctuations induced by changing the source strength itself. Therefore using the model of classical emissions 
^,[H||| we should be able to predict the spectrum of rapidity fluctuations at particle formation time. 

The nature of the fluctuations in the source density is complicated by the fact that the source density 
itself satisfies a renormalization group equation |^],[|, and is correlated with what the source would have 
been at other values of rapidity. One of the purposes of this paper is to disentangle this dependence and 
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to predict the general from of the multiplicity fluctuation spectrum. In Sect. II we shall show that as a 
consequence of the renormalization group behavior, the fluctuation spectrum exhibits the Koba, Nielsen, 
and Olesen (KNO) scaling [^8|. The KNO scaling for hadron-hadron collisions states that the probability 
of producing a given multiplicity of particles N at some given high energy E, which we denote dP/dN, 
multiplied by the average multiplicity of produced particles at this energy N(E) is described at all energies 
by the same scaling function of N/N(E) ^| 

N(E) ^ = f(N/N(E)), (6) 

characteristic of given hadron species. That way the energy dependence of the high energy inclusive cross 
sections comes in only through the average multiplicity of the produced particles N(E). 

The nature of this scaling is slightly modified in the case of a nucleus since in the nuclear collisions 
there is another parameter in the problem — atomic number of the nucleus A. Many nucleons in the 
nucleus enhance multiple interactions, changing the shape of the KNO distribution. In this paper we will 
for simplicity consider only the collisions of identical nuclei having the same atomic number A. The relative 
width of the distribution depends on the baryon number of the nucleus, and must scale roughly like 1/V A s , 
where S could be 2/3, 1 or 4/3 depending on the saturation model, as will be discussed below in Sect. II. 
This scaling arises from the fact that nucleons separated by transverse distances greater than or of the order 
of a Fermi from one another must act independently. Thus for nuclear collisions we predict that the KNO 
scaling function will depend on A 

N(E)^ = f(A,N/N(E)). (7) 

We will derive this result in Sect. II. It implies that for a nucleus of fixed size A one should observe KNO 
scaling, but the scaling function will differ from a nucleus to a nucleus. We will also argue that the form of 
this distribution is 



/ ~ exp 



-const A s {\Jn/N{E)~ if 



(8) 



[see Eq. (pl|)]. This multiplicity distribution has a width 

SN ~ N(E)/Vas. (9) 

Of course the width SN is parametrically of the order of y/A?, since, as will be shown in Sect. II, N(E) ~ A s , 
but the signal for the KNO form is its dependence on energy through the total multiplicity. 

Of course, due to the final state interaction of particles, the number of particles in the final state is 
not expected to be much changed from the initial conditions. We expect therefore that fluctuations in the 
total particle multiplicity as a function of rapidity on an event by event basis at the particle formation time 
should be reflected in the final st ate dis tribution of produced particles. This should be true up to at least 
some number such as a s times y dN/dy. If we look for very rare fluctuations, then such fluctuations in the 
final state interactions should not obscure a huge initial state fluctuation. 

In Sect. Ill we will derive another feature of the underlying classical dynamics which manifests itself in 
the two particle distribution function. To compute the two particle distribution, for each event we take the 
multiplicity in some bin at y\ and multiply by the multiplicity at another rapidity bin y 2 . If the classical 
fields are the ones responsible for the distribution of the produced particles, we shall show that at the leading 
order in A the two particle distribution function factorizes into 

, dN dN 

D(y u y 2 ) = — —. 10 
dyi dy 2 

The factorization is demonstrated in Fig. [| Of course this factorization would happen if the distributions 
were uncorrelated at rapidities y\ and y 2 . We shall argue that they are in fact very tightly correlated in 
Sect. IV by introducing a correlation function which would allow one to disentangle the classical effects from 
other possible mechanisms of particle production. 
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The classical field responsible for this distribution is constant over a large rapidity interval. Therefore if 
we measure the multiplicity fluctuation in a rapidity interval centered around y\ which is several sigma away 
from the average, the multiplicity in a neighboring rapidity interval centered around 2/2 should be roughly 
the same multiple sigma fluctuation from the average. This very remarkable correlation is a measure of the 
classical coherence of the field, and is a direct measure of the underlying strong like dynamics of the gluon 
field. We shall make this argument firm in Sect. IV by introducing a correlation coefficient at a fixed impact 
parameter of the nuclei C(b, 2/1,3/2) in Eq. (32). The coefficient C (6, 2/1,2/2) is an experimentally measurable 
observable which is equal to 1 when the particle productions at the rapidities 2/1 and 2/2 are correlated on the 
event-by-event basis, and is less than one otherwise. We will argue that if the classical picture of emissions 
is true this correlation coefficient should be 1 over large intervals in rapidity, 1 2/1 — 2/2 1 ~ 1/ck s . It should fall 
off for wider rapidity intervals. This prediction could be checked experimentally at RHIC and LHC. 

In Sect. IV we also discuss the differences and similarities of the results of our classical approach and 
the conclusions for rapidity correlations which could be driven out of the old reggeon theory. 

We summarize our results in Sect. V. 



II. NATURE OF KNO SCALING 



In this section, we will study the fluctuation spectra of produced particles. Strictly speaking, we are 
studying the fluctuation spectrum of produced gluons in early stages of the collisions. The number of gluons 
is of course modified by subsequent interactions with other gluons and in their subsequent transmutation 
into pions. It is expected nevertheless that the number of produced pions is close to the number of initially 
produced gluons. This is because the gluons are thermalized largely by two body collisions in the early 
stage of the collision, since the coupling is weak if the typical transverse momentum of the gluons is large. 
This is expected if the multiplicity per unit area A 2 of Eq. (Q) is large compared to Aq CD as should 
be the case for large nuclei at asymptotically high energy. The two body collisions change the transverse 
momentum distribution of gluons, but preserve the total gluon number. In this paper we are interested in 
the multiplicities of gluons integrated over all transverse momenta (or coordinates), which do not change 
under thermalizing two gluon collisions. At later times after thermalization, we expect that entropy will be 
approximately conserved. At late stages, the entropy is converted into pion number as the system cools. 
We expect therefore that dN g i uon /dy ~ dN P i on /dy. In general there may be some weak dependence of this 
proportionality upon multiplicity, but for our purposes this weak dependence will not be important. 

Therefore, although we compute the initial multiplicity fluctuations in gluons, this should be reflected 
in the multiplicity fluctuations of produced pions. 

In the McLerran-Venugopalan model of the small x hadronic wave function, one computes a classical 
field which arises from a source density J^j5|. This source density is at rapidities much larger than that of 
the classical field which we compute. Eventually the source density is itself integrated over to be included 
in the source for the next series of classical emissions |||| , as was described in the introduction. 

Let us first compute the fluctuation spectrum of produced particles for a fixed source. In this case, the 
wavefunction corresponding to this classical field is a coherent state. Coherent state wavefunction generate 
Poisson distribution in multiplicities. Therefore the typical fluctuation scale in a unit of rapidity is of order 
y/dN /dy. We will soon see that this scale is small compared to that generated by the fluctuations in the 
source itself. 

To understand the fluctuations induced by the source of the classical fields, we construct a model which 
has most of the features of BFKL evolution and fluctuations in the source density. We use intuition developed 
from understanding the renormalization group structure of the McLerran-Venugopalan model. We use the 
equation 

^ = kN( V ) + v^iV(^C(2/)- (11) 
dy 

In this equation, N(y) is the total number of particles in the rapidity range between ytarget and 2/, 

f v dN 
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The first term on the right hand side of the evolution equation ( |ll| ) is the toy model analog of the kernel of 



the BFKL evolution equation (thus k = ctp — 1 



4a. JV, In 2 



is the intercept of the BFKL pomeron M ) . For 



a series of classical emissions one could conclude that the number of produced particles is proportional to 
the total number of partons off which the new particles are emitted. This is reflected in the first terms of 
Eq. ( |Ti"| ) . We note that due to the classical emission picture we can write Eq. (|lTJ) for particle multiplicities 
and not just for the cross sections (The BFKL equation was originally written for the cross section.). The 
reason for that stems from the fact that at the high energies considered here the total cross sections become 
independent of energies, thus making the energy dependence of inclusive cross sections and multiplicities 
identical. If we ignore the second term on the right hand side of Eq. ([ll]) the solution of the equation would 
be 



ay 



(13) 



which is our analog of the solution of the BFKL equation. Comparing Eq. (13) to the usual one BFKL 
pomeron exchange result || we note here that No ~ a 2 . The A-dependence of iVo is determined from the 
predictions for multiplicity inside of the saturation region. All saturation models [p UT^ , |29[| agree that the 
multiplicity of the produced particles in the saturation region is given by 



ay 



(14) 



with Q s (y) the saturation momentum, which is denoted by A in Eq. (||). The factor of A 2 I Z in Eq. ( |l4| ) 
results from the integration over the impact parameters of the nuclei. In the Glauber-Mueller-type saturation 
models f30|| the saturation momentum scales as Q 2 ~ A 1 / 3 due to multiple rescatterings within the nucleus. 
In the approaches based on the BFKL equation ||,|2^,^l[ Q 2 ~ A 2 / 3 . In the Reggeon approach [pi|-|26|] Q 2 
is independent of A. We can summarize all these results by writing Nq ~ A s , with 6 being dependent on 
the particular model at hand. In Glauber-Mueller model S = 1, in BFKL-based approaches 8 — 4/3 and in 
Reggeon calculus 5 = 2/3. 

The second term on the right hand side of Eq. ([ll]) is stochastic and in our model describes the 
fluctuations induced at each stage of the evolution equation. The probability distribution of the fluctuation 
of the stochastic term £ is Gaussian 



z = J[d(} e -^f dyc2(v) . 



(15) 



Its origin is from the renormalization group || . At each step a classical field is induced, which when small 
fluctuations are computed gets converted into a source for the next step in the evolution. Of course the 
classical field itself has Poissonian fluctuations which near the center of the distribution should be Gaussian. 
Therefore, the induces source at the next step will have fluctuations built in which are correlated with the 
fluctuations in the previous step. If we look over an interval of unit width in rapidity, the wei ght of these 
fluctuations is of order dN/dy. Our stochastic source C is weighted by the function Z of Eq. ( |l5| ) so that 
its fluctuations in one unit of rapidity are also of order one. Using the fact that k N ~ dN/dy which is true 
for the solution ( |l3| ) of Eq. ( jll| ) without the stochastic term, and, as we shall see below, is approximately 
valid for the full Eq. ([ll]) , we can replace the weight of the stochastic fluctuation term by \J dN/dy — > V V N 
arriving at the evolution equation (11). We also note that k! ~ a s . 
The solution to Eq. ([ll]) for fixed ( is 



N(y) 



N(y) 



InEq. (M 



/ . ' / dy e 

o 



Civ') 



N(y) = e KV N Q 



(16) 



(17) 



where we will in the future set ytarget = 0. Note that N is the typical average multiplicity at fixed y which 
one would have in the limit when fluctuations are turned off. 



G 



We can now compute the multiplicity distribution function dP/dN . To this we must integrate over the 
fluctuations fields C with the constraint that the total multiplicity is given by Eq. (|f6|). We have 



§f = l / K]6XP 



-5 / dy'( 2 (y') 



< s I N(y) - [y N(y) - -y- 



dy e 



(18) 



In this integral, we replaced the upper limit of integration in the expression for N(y) by infinity, since we 
are typically interested in large values of y. This makes the analysis simpler in what follows. 
To evaluate the path integral, one should decompose 



C(y) 



e Ky/2 c n L n (Ky), 



(19) 



where the L n 's are Laguerre polynomials lP n . Using the orthogonality condition for Laguerre polynomials 
one can see that all the integrals over c n 's in Eq. (|l^) become Gaussian with the exception of the integral 
over Co in the numerator, which is fixed by the delta function. The integrals over c„'s for n > 1 can be done 
in closed form and are canceled by the same integrals in Z in the denominator. The final answer is given by 
the Co integration. 

In the end, we find that 



itf/ \ dP 

N(y) 

y!JJ dN 



KN N(y) 
2-kk'N 



exp 



kN 



The multiplicity of the produced particles Nq is very large. 
N = ~N. We obtain 



1 - 



/ N 
N(y) 



(20) 



That allows us to expand Eq. ( BO) around 



— , x dP 
~dN 



KNn 



2ttk 



7 exp 



2k' 



N 

1 - — 



N(y) 



(21) 



This result for the multiplicity distribution is simple to understand. The evolution equation connects the 
produced particle to some initial spectrum of fluctuations. The relative importance of these fluctuations de- 
creases as the multiplicity increases. Therefore the fluctuations within an interval of width l/a s in rapidity 
set up the fluctuations at higher values. This is the factor of kNq/k' in the exponent of Eq. (|2l|). Remem- 
bering that k ~ k' ~ a s and N ~ a 2 s A 6 we conclude that the relative width of the multiplicity distribution 
scales like l/a s \/~A^. The fluctuations are all correlated all the way down the chain and therefore should 
only depend upon the ratio N/N(y). This dependence is the essence of KNO scaling. The only important 
variable is the multiplicity divided by the average multiplicity. Of course there is also a dependence on Nq 
because there are more independent emitters for a nucleus than for a hadron. This is reflected in the fact 
that N - A 5 . 

Note that the width of this distribution is 



SN 2 ~ N {y)K /N n . 



(22) 



Parametrically <5iV 2 is linear in A 8 , as N — A 5 and N — N — A s . 

In general, the fluctuation spectrum for the first few emitters which generate the fluctuations in the 
distribution may be different from Eq. ( |2lh. The exact shape of the KNO scaling function will probably 
differ from the one given by Eqs. ( |20| ) and (|2l|). Nevertheless, the physical picture we have generated still is 
true, and therefore we expect that in general the distribution will be of the form 



N(y)^ = f(N ,N/N(y)) 



(23) 



and that the typical width of the distribution will scale with A and a s in the same way as the width in our 
model given by Eq. (|22| ) . The physical picture we have of the small x gluon distribution functions therefore 
automatically has KNO scaling [cf. p3|]. 



7 



In terms of the strong coupling constant the relative width 5N 2 / N 2 {y) ~ l/ct^ (see Eq. (p2|)). It has 
been argued in the framework of McLerran-Venugopalan model that in high energy nuclear scatterings the 
high parton density sets the scale of the running coupling constant If this assumption is true, than 

one could conclude that at very high energies No ~ a 2 {A ta rget), since Nq is the multiplicity of the particles 
in the fragmentation region. A target is the transverse momentum scale characterizing the target nucleus. At 
the same time k ~ a s (A), as it determines the evolution of the multiplicity N(y) at some large rapidity y. 
(A is given by Eq. (||).) k' depends on y through its dependence on A~ N(y). However, if we allow k' to 
depend on y then, since the integral over y' in Eq. ( |l6| ) is dominated by y' ~ one can see that the most 
important value of k' would be at the fragmentation region. Thus k' ~ o: s {^-target)- In the light of the above 
estimates and using Eq. (|22] ) we conclude that the relative width depends on the average multiplicity 

6N 2 1 1 

^ ^ (24) 

N 2 (y) "s(A) a s {K target ) a s (N)' 

The high multiplicity of the produced particles would create a large momentum scale A, which would make 
the coupling constant a s (A) small. Thus the relative width of the KNO distribution given by Eq. (p4| ) will 
get larger. Moreover, since the change in width would depend on the average multiplicity of the produced 
particles that would violate the KNO scaling. Thus our prediction of KNO scaling will start to break down 
slowly at very high energies. At the same time since the rate of change of the running coupling constant 
a s (A) slows down at large A the violation of KNO scaling will decrease as the energy gets higher. 



III. CLASSICAL NATURE OF THE DISTRIBUTION FUNCTION 



Another feature of the classical field particle production arises from considering the distribution function 
D(yi, 2/2)- We define this to be the two particle distribution function measured in the following way: In each 
event or class of events, measure the multiplicity of particles at rapidities y± and y 2 in bins of width dy\ and 
dy 2 . Multiply these multiplicities together to generate 

d 2 N(y 1 ,y 2 ) = D(y 1 ,y 2 ) d Vl dy 2 . (25) 

Then average D(yi, y 2 ) over events or classes of events. 



sources of nucleus #1 




A class of events might be generated by putting a cut on say the total multiplicity of particles in the 
neighborhood of zero rapidity. In the class of events where this cut is satisfied, one could measure the 
multiplicity in different regions of rapidity and average over the class of events. 

The distribution function D may be generated as an expectation value of the number operator 
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D{yi,y 2 ) = J {a\yi,PT)a{yi,p T )a){y2,qT)a(y 2 ,qT)) ■ (26) 

In the classical field limit, this expression is of the form 

(a\y 1 ,p T )a(y 1 ,p T )a\ y2: q T )a(y 2 ,q T )) ~ A i {jp)A i {-p)A i {q)A i {-q) (27) 



where A 1 is the classical field produce by the color sources in the colliding nuclei |4|,|9[|1C|] . When one averages 
over the sources which generate the fields, one finds that 

, /dN\ /dN\ 

DM "(*)(*;)■ (28) 

This is true up to corrections which are of order A~ 2 / 3 , that is of the order of one over the area of the nuclei. 
This statement could be understood from Fig. |3[ The classical fields producing the gluon connect to several 
color sources. In the leadingorder in A the fields generating gluons at y\ and y 2 connect to different sets of 
color sources, making Eq. ( pq) true. 

Note that this tells us that the fields are essentially trivially correlated in the longitudinal direction. 
The connected piece of D(yi,y 2 ) has vanished entirely, up to corrections which go like one over the area of 
the nuclei. This lack of correlation is not by itself evidence of a classical field. It could occur if there were 
for example no correlations at all in the longitudinal space. The structure of the KNO distribution itself 
suggests that this is not the case. In the next section, we shall see that there is a correlation function which 
dramatically illustrates the classical correlation. 

IV. CLASSICAL CORRELATION 

If there are classical fields, then, as we have seen above, due to boost invariance, these classical fields are 
independent of rapidity over a wide range of rapidity [QJ^] . Therefore the multiplicity should on the average 
be the same over the range of rapidity where the classical field theory is valid. This is typically of order l/a s 
in the classical field approach. At the same time a s should be small when the density of produced particles 
is large making this correlation length in rapidity large. 

This effect can be measured in the following way: Measure the rapidity density in some bin of width dy\ 
around y\ and dy 2 around y 2 . Require that dyi and dy 2 are large enough so that the statistical fluctuations 
in the rapidity in a given event are small compared to the total multiplicity. From our picture of classical 
particle production it follows that the multiplicity around y\ should be the same as around y 2 up to statistical 
fluctuations. 

If we applied this analysis to an average event, the result above would be trivial. On the other hand, if 
we look at an event where the fluctuation is very rare in the bin around t/i, we predict that there will also be 
the same rare fluctuations around y 2 \ Rare events tend to fluctuate together over a wide range of rapidity! 

A measurement of the width of the correlation length tells us something about the underlying classical 
dynamics and is interesting in itself. 

First we note that in a nuclear collision 2L. is a function of rapidity y, impact parameter of the nuclei b, 
and of the configuration of color charges in the nuclei in this particular collision, which we will symbolically 
denote p: 

dN dN , , , 

— = —(y,b,p). (29) 

ay dy 

If one measures in a number of events and then takes the average value at a fixed impact parameter of 
the colliding nuclei b the result should correspond to the averaging of the theoretical prediction for ^ over 
all configurations of color charges in the colliding nuclei p, which we write as 

f <*M) ■ (30) 



9 



Event by event fluctuations in ^(y, b, p) are characterized by the variance of that quantity which we will 
define as 



V 



dN] (dN, , A 2 \ /dN, , A 2 

*\ = \W vAp) ))-W vM ) r (31) 



Thus we can also introduce the correlation function of the numbers of particles measured at two different 
rapidity points yi and y 2 in several events with the same impact parameter (that is for several configurations 
of color charges p). According to the standard mathematical methods we define the correlation coefficient 

C(b;y 1 ,y 2 ) = , = • (32) 



V 



dN 



din 



V 



dN 



dyi 



We can neglect the y dependence of ^(y, b, p) for the rapidity intervals of the size 1/a. That way one can 
see that when \yi — 2/2|<l/a we can neglect the variation of ^{b,p) with y, which is purely statistical 
in that interval, and put rj J^-, which leads to C(b; yi, y 2 ) = 1- This is the prediction of the classical 
emission picture discussed above and in the correlation function C(b; 2/1,2/2) should be close to one for 

\yi ~ 2/2 1 l/ct. When \yi — t/2|sJl/o; the y dependence in in the interval between yi and 2/2 becomes 
important and and become uncorrelated, which would lead to C(b; yi, y 2 ) being smaller than one. 

One might note that in the leading powers of A approximation both numerator and denominator in 
the expression for C(6; 2/1,2/2) given by Eq. (^) are zero. The non-zero terms appear once we include the 
corrections which are subleading in A. The first non-zero terms in the numerator and denominator of Eq. 

dN 



([32|) are suppressed by A 2 / 3 compared to the leading contribution in, for instance, , (b,p). That is 



dyxdyi 

if we calculate d y^ y2 (b, p) in the classical approximation we would obtain the leading term in A given by 
Eq. (p8|), which is canceled in the numerator of Eq. (|3^), and the subleading terms, which are suppressed 
at least by A 2 / 3 , but still arise from the classical approximation. The largest of these subleading terms is 
suppressed by exactly A 2 ! 3 and results from the situation when both particles at y\ and 2/2 are emitted from 
the same nucleon. Even though it is subleading in A this term is nevertheless larger than the correlations 
induce by quantum corrections connecting the fields that produce the particles at y\ and 2/2 , since the latter 
are also suppressed by extra powers of a s - The quantum correction within one of the fields do not violate 
the factorization of Eq. ( p8|) and, therefore, cancel in the numerator of Eq. (^). Thus the fact that the 
correlator in Eq. (j3^) is equal to 1 because of the terms in the numerator and denominator which are 
subleading in A does not influence the fact that the correlations are classical. 

In order to understand Eq. ([32]) better let us consider different possible definitions of the correlation 
function. For rapidity correlation function defined as 



/^{b,p) -fL(b,p)fL( bjP ) 

-rp I.. „, \ _ dyidy-2 \ i r I dyi V ' ^' dyi V ' rl , . 

KM,9 ' ) - fiMfH 



we expect that 



n p ( yi ,y 2 ) = 0{A—z . (34) 



Note, that in Eq.(|32|) we fix the value of the impact parameter (transverse distance between centers of 
colliding nuclei) which is also a characteristics of the event. The A~ 3 - corrections stem from region of 
integration px ~ <Zt in Eq. (pjj). This is the subleading correction to the factorization of Eq. (p8|). It could 
also be viewed as resulting from the contribution when the two fields of Fig. || share the same source of 
color charge. That term would definitely be suppressed by A 2 ! 3 . Thus the correlation coefficient of Eq. ( |33| ) 
defined without averaging over events does not reflect the correlations that we are interested in. However, R p 
is not the correlation coefficient usually employed in the experiments. One defines the correlation function 

"^average by 
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1 / da \ 1_ / da_ \ 1 / da_ 

a tot \d yi dy 2 / <T tBt \dyi/ b <T t o t \dy 2 , v 
Average (2/1, 2/2) = 7 \ J—\ ( 35 ) 



1 / da \ 1 / do - 



^tot \dvi/ biP "^t \dy2, bp 

where <Jtot, jffj and dy f^ V2 are total, single and double inclusive cross sections for nuclear collisions that 
are defined by averaging over all events (p) and integrating over all impact parameters (6*). As was shown 
long ago in the Pomeron approach function lZ aver age generally is of the order of unity even in the case if 
have IZp = [fl6|| . However, we would like to point out that the main correlation which makes lZ average 
large is the simple correlation in impact parameters which has a very simple underlying physics, which 
states that partons (gluons) are produced independently but at the same value of impact parameter. In 
McLerran-Venugopalan model one might expect that at high parton densities or better to say in collisions 
of heavy nuclei the b t - distribution for all physical observables are the same, namely, Q(Ra — bt). In this 
case IZaverage — > showing a classical field emission even without the event selection (b t fixing) . If one does 
not integrate over the impact parameter in Eq. ( |35| ) then without the impact parameter correlations this 
correlation function will be small again failing to reflect the type of correlations we are interested in. 

The general picture of correlations for high parton density QCD turns out to be very similar to pattern of 
correlations predicted by the Reggeon approach which has been discussed in details for three decades 27] . 
The first observation is that in both approaches Eqs. (^7|) and ( |33| ) lead to conclusion that the secondary 
gluons (hadrons) are originated from the independent production of clusters with the mean multiplicity 
(N) = J ^{y,p,b t ) dy. Therefore, the probability to emit k x (N) gluons or, in other words, k clusters, is 
equal to 

Pk = e -W . (36) 



Note, that Eq. (36) is a well known Poisson distribution. 



The second observation stems from the classical field emission, namely ^jp(y 5 P,bt) — d — Const(y) 3> 1 
for the rapidity interval of the order of l/a s 3> 1. We will discuss why d> 1 below. Assuming this we can 
easily see that we can neglect statistic fluctuations in N(y, p,b t ), and, therefore, (N) — Yd where Y is the 
accessible rapidity interval (Y oc \/ot s ). We can claim even that N(y,p,bt) which is defined as the number 
of gluons ( hadrons) in the rapidity interval y — Ay y + Ay is equal to 

N(y,p,b t ) = kx d x 2Ay . (37) 

Therefore, measuring N(y, p, b t ) we fix the number of clusters k or, in other words, we select a configuration 
which is produced with probability (see Eq. (|36])). 

Let us recall the main predictions that will hold both in classical field approach ( McLerran-Venugopalan 
model ) and in the Reggeon description. If we select the events where the multiplicity of the particles in 
the rapidity bin y\ — Ay -j- t/i + Ay around yi (N(yi)) is fixed and average the multiplicity in the rapidity 
interval y2 — Ay y2 + Ay over these events we predict that it will be equal to 

(N(y2)) N(yi)flxed = N( yi ) . (38) 

Also the number of neutral pions in rapidity interval y\ — Ay -j- yi + Ay (N (yi)) averaged over the events 
with a fixed number of charged pions in the same bin should be equal to 

(N (yi)) N * Hyi)fixed = \ N ch (y x ) , (39) 

where N ch (yi) is the number of charged pions. Both models predict that for heavy ion- ion collisions in the 
kinematic region of saturation PJl5[| we have 

TlaverageiV 1,V2) > (^J^j ■ (40) 

One can see that all these prediction are direct consequences of Eq. (|3^) and Eq. ( |37"|) . Much more detailed 
predictions could be found in Refs. pTVpTl]. 



11 



The high density QCD and soft Reggeon approaches can be distinguished by measuring the transverse 
momenta distributions or/and transverse momenta correlation between produced particles. Indeed, in high 
parton density QCD the typical transverse momentum of produced gluons (saturation momentum Q s (A, x) ) 
is large and depends on A PJD| while in the Reggeon approach this momentum is constant and rather small 
- about 2 GeV, which results from the slope estimates of the soft pomeron. However, at first sight, this 
main difference does not influence the shape of the rapidity correlations and only increases the multiplicity 
of produced particles making our predictions more reliable in this case. Indeed, the saturation of the gluon 
density leads to d oc Q 2 S (W : A) M^P^f! where W is the energy in the center of mass frame and 
Ql(W, A) is the saturation scale Eq. (|14|), which is similar to A of Eq. (||). Since Q 2 (W, A) increases with W 
and A (at least Q 2 S oc A* |,|,|Jl|j2| but it could even be proportional to is |3lJ] ) we expect d> 1 which 
makes all of our estimates much more accurate than the results of the Reggeon approach 122 
as the shape of the rapidity correlations is concerned the saturation of the parton densities 
leads to sufficiently large correlation length which is proportional to l/a s (Q s ). For heavy nuclei and/or 
high energies Q s increases and a s — > 0. This fact leads to Eq. (|3^) being better justified in the high parton 
density QCD than in the case of the Reggeon approach. 




V. CONCLUSIONS 



In this paper we have considered the consequences of the model of classical emissions on the particle 
production mechanisms in heavy ion collisions ||. We have argued that the classical fields do not vary 
significantly over the rapidity intervals of the order of l/a s . This led us to conclude that the multiplicities of 
the particles produced in the early stages of the collisions (gluons) are correlated over l/a s units of rapidity 
on the event-by-event basis. We then gave an argument demonstrating that the number of gluons generated 
in a particular collision is proportional to the number of pions produced in the final state. Therefore the 
correlation in the multiplicities of the produced gluons would reflect themselves in the multiplicities of the 
final state pions. We have then constructed the correlation function C(b;yi, 2/2) in Eq. (|3|), which could be 
measured experimentally and is equal to 1 when the multiplicities of particles at the rapidities 2/1 and 2/2 
are correlated and is less than 1 otherwise. Surprisingly this effect comes from the terms in the numerator 
and denominator of the expression for C(b; 2/1,2/2) given by Eq. (j32^ ) which are subleading in A. We predict 
C(b; y±, 7/2) to be close to 1 when \yi — y<i\ < l/a s and fall off for larger rapidity intervals. 

We have also used the model of classical particle production to explain KNO scaling. In a simple toy 
model which resulted from trying to mimic the main features of the classical emissions we have reproduced 
KNO scaling for the multiplicities of produced particles. The result is given in Eq. (|2l]). Even though 
the exact shape of the KNO distribution function is, probably, different from our toy model prediction, we 
believe that the model captures its main features. For the case of nuclei we predict the KNO function to 
depend on the atomic number A in addition to the usual dependence on N/N(y). We predict that at very 
high energies KNO scaling in nuclear collisions will be violated due to the running of the coupling constant, 
which would get smaller as parton density increases. 

The main results of this paper can be summarized in the following way: 



1. We have derived KNO scaling from the classical emission picture of particle production (see Eqs. ( ppj ) 
and (§!])). 

2. We predict rapidity correlations on the scales Sy ~ l/a s in the particle (pion) production in heavy ion 
collisions on the event-by-event basis (Sect. IV). 



3. We proposed a correlation coefficient C(&; 2/1,2/2) which would allow one to measure the predicted 
correlations (see Eq. (fj2|)). 
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